#pragma once

#include "matrix3.h"
#define  TEMPLATE template<typename T>
#define  MATRIX3 Matrix3T<T>


BEGIN_NAMESPACE_OXGE


TEMPLATE
const MATRIX3 MATRIX3::IDENTITY(1,0,0,
								0,1,0,
								0,0,1);

TEMPLATE
const MATRIX3 MATRIX3::ZERO( 0,0,0,
							 0,0,0,
							 0,0,0 );



//m * m
TEMPLATE
MATRIX3 MATRIX3::operator * (const MATRIX3& matR)const
{
	return MATRIX3( RowVector(0).dot(matR.ColVector(0)), RowVector(0).dot(matR.ColVector(1)), RowVector(0).dot(matR.ColVector(2)), 
					RowVector(1).dot(matR.ColVector(0)), RowVector(1).dot(matR.ColVector(1)), RowVector(1).dot(matR.ColVector(2)), 
					RowVector(2).dot(matR.ColVector(0)), RowVector(2).dot(matR.ColVector(1)), RowVector(2).dot(matR.ColVector(2)) ); 
}
//determination
TEMPLATE
T MATRIX3::Det()const
{
	return m00 * GetSubMatrix(0,0).Det() - m01 * GetSubMatrix(0,1).Det() + m02 * GetSubMatrix(0,2).Det();
}

//m ^ p
TEMPLATE
inline MATRIX3 MATRIX3::operator^(int p)const
{
	MatrixType ret = MatrixType::IDENTITY;
	for(int i=0;i< abs(p) ; ++i)
		ret *= (*this);

	if( p < 0 )
		return ret.GetInverse();

	return ret;
}

//sub matrix
TEMPLATE
typename MATRIX3::SubMatrixType MATRIX3::GetSubMatrix(int r, int c)const
{		
	T matSub[SubMatrixType::NCol * SubMatrixType::NRow];		
	size_t i2 = 0;
	for(size_t i=0;i< NRow; ++i)
	{
		if( i==r )continue;
		for(size_t j=0;j<NCol;++j)
		{
			if( j==c )continue;
			matSub[i2] = (*this)(i,j);
			++i2;
		}
	}		
	return SubMatrixType(matSub); 
}


TEMPLATE
bool inverse(MATRIX3& dest, const MATRIX3& src)
{
	T detsrc = src.Det();
	if( detsrc == 0 )
	{
		return false;
	}

	T rDet = T(1)/detsrc;

	if( &dest != &src )//dest not the src
	{
		for(size_t i=0;i<src.NRow;++i)
		{
			for(size_t j=0; j<src.NCol;++j)
			{
				dest(j,i) = ((i+j)%2==0 ?(rDet):(-rDet)) * src.GetSubMatrix(i,j).Det() ; 
			}
		}
	}
	else
	{
		Matrix3T<T> srccpy = src;
		for(size_t i=0;i<srccpy.NRow;++i)
		{
			for(size_t j=0; j<srccpy.NCol;++j)
			{
				dest(j,i) = ((i+j)%2==0 ?(rDet):(-rDet)) * srccpy.GetSubMatrix(i,j).Det() ; 
			}
		}
	}
	return true;
}

//global I/O
TEMPLATE
std::ostream& operator<<(std::ostream& os, MATRIX3 mat)
{
	size_t iRow = 0;	
	for(;iRow < mat.NRow-1;++iRow){
		size_t iCol = 0;
		for(;iCol < mat.NCol-1 ; ++iCol){
			os << mat(iRow, iCol) << "\t";
		}
		os << mat(iRow, iCol);
		os << '\n';
	}
	size_t iCol = 0;
	for(;iCol < mat.NCol-1 ; ++iCol){
		os << mat(iRow, iCol) << "\t";
	}
	os << mat(iRow, iCol);
	os << ';';

	return os;
}


//unit test
TEMPLATE
void MATRIX3::UnitTest()
{
	MATRIX3 mat;
	MATRIX3 mat1 = ZERO;
	mat1 = IDENTITY;
	mat.Det();

	mat^3;
	mat * mat1;

}


END_NAMESPACE_OXGE


#undef TEMPLATE
#undef MATRIX3
